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Abstract: Combustion instability modeling of 
Solid Rocket Motors (SRM) remains a topic of 
active research. Many rockets display violent 
fluctuations in pressure, velocity, and temperature 
originating from the complex interactions 
between the combustion process, acoustics, and 
steady-state gas dynamics. Recent advances in 
defining the energy transport of disturbances 
within steady flow-fields have been applied by 
combustion stability modelers to improve the 
analysis framework [1, 2, 3]. Employing this more 
accurate global energy balance requires a higher 
fidelity model of the SRM flow-field and acoustic 
mode shapes. The current industry standard 
analysis tool utilizes a one dimensional analysis 
of the time dependent fluid dynamics along with 
a quasi-three dimensional propellant grain 
regression model to determine the SRM ballistics. 
The code then couples with another application 
that calculates the eigenvalues of the one 
dimensional homogenous wave equation. The 
mean flow parameters and acoustic normal modes 
are coupled to evaluate the stability theory 
developed and popularized by Culick [4, 5]. The 
assumption of a linear, non-dissipative wave in a 
quiescent fluid remains valid while acoustic 
amplitudes are small and local gas velocities stay 
below Mach 0.2. The current study employs the 
COMSOL multiphysics finite element framework 
to model the steady flow-field parameters and 
acoustic normal modes of a generic SRM. The 
study requires one way coupling of the CFD High 
Mach Number Flow (HMNF) and mathematics 
module. The HMNF module evaluates the gas 
flow inside of a SRM using St. Robert’s law to 
model the solid propellant burn rate, no slip 
boundary conditions, and the hybrid outflow 
condition. Results from the HMNF model are 
verified by comparing the pertinent ballistics 
parameters with the industry standard code 
outputs (i.e. pressure drop, thrust, ect.). These 
results are then used by the coefficient form of the 
mathematics module to determine the complex 
eigenvalues of the Acoustic Velocity Potential 
Equation (AVPE). The mathematics model is 
truncated at the nozzle sonic line, where a zero 


flux boundary condition is self- satisfying. The 
remaining boundaries are modeled with a zero 
flux boundary condition, assuming zero acoustic 
absorption on all surfaces. The results of the 
steady- state CFD and AVPE analyses are used to 
calculate the linear acoustic growth rate as is 
defined by Flandro and Jacob [2, 3]. In order to 
verify the process implemented within COMSOL 
we first employ the Culick theory and compare 
the results with the industry standard. After the 
process is verified, the Flandro/Jacob energy 
balance theory is employed and results displayed. 

Keywords: Combustion Instability, Solid Rocket 
Motor. 

1. Introduction 

Combustion instability modeling within the 
Solid Rocket Motor (SRM) industry has been 
largely limited to the JANNAF Standard Stability 
Prediction code (SSP) embedded within the Solid 
Propellant Performance program (SPP’04) [5]. 
The SPP’04 software is the industry standard for 
modeling SRM internal ballistics and contains 
and extensive, and ever growing, array of 
modules for predicting various mechanisms 
including erosive burning, boundary layer loss, 
grain regression, etc. The SSP code, originally 
developed solely for axial instabilities, is based on 
the Culick [4] acoustic formulation for 
combustion instability, and utilizes internal 
ballistics parameters derived by the SPP’04 
analysis. Recent advances in combustion stability 
modeling by Flandro et al. [6] employ an energy 
based approach, allowing for the inclusion of 
rotational and thermal oscillations. 
Independently, Myers [1] developed a general 
unsteady energy corollary for the transport of 
disturbances in steady flows. Application of the 
Myers corollary to combustion instability 
modeling was developed and derived in detail by 
Jacob [3]. With the increase in computational 
resources and a more detailed theory for stability 
modeling the utilization of a high fidelity analysis 
tool, like COMSOL mutiphysics, is now practical 



and tractable for the SRM industry. The present 
study investigates the feasibility and advantage of 
employing COMSOL in the prediction of 
combustion instability in SRMs. 

1.1 Acoustic Wave Equation Stability Model 

Beginning with a linearized form of the 
governing equations, a wave equation can be 
constructed in a manner similar to classical 
acoustic theory. The assumption is made that the 
steady portion of the flow parameters do not 
change over the time period of a single acoustic 
oscillation. This allows for the flow parameters to 
be written as a sum of mean and oscillating 
parts, Q = Q + q\ respectively. After neglecting 
the terms of second order in the oscillatory 
parameter, the linear inhomogeneous wave 
equation becomes, 
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The traditional approach to calculating linear 
stability starts by first assuming the classical 
acoustic formulation, h = 0, and solving for the 
normal modes. The difference of the homogenous 
and inhomogeneous solution are compared 
through a spatial average. Further detail can be 
found in the works of Culick [4]. Assuming that 
the oscillatory parameter can be modeled as a 
damped oscillator, 
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the shift in frequency, co n> and growth constant, 
a n , can be determined. An alpha greater than zero 
represents a linearly unstable system and a growth 
constant less than zero depicts a stable system. 
The equation for alpha represents a list of sources 
and sinks of acoustic energy within the 
combustion chamber. Written in the SSP 
nomenclature, the main contributors are, 
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1.2 Energy Corollary Model 


Investigating combustion instability through 
the lens of an energy equation allows for a more 
complete accounting of all flow disturbances. The 
historical focus on acoustic oscillations due to the 
often audible tones from unstable combustors 
captures a majority of flow unsteadiness, but miss 
contributions made by vorticity and entropy 
oscillations. These contributions, which may be 
small in some configurations, become important 
when unsteady combustion and complex 
geometries are considered. Myers [1] defines an 
energy corollary by employing the continuity, 
momentum and entropy equations and expanding 
the flow variables into slow changing and fast 
changing time scales. 
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Jacob [3] applied the Myer framework to SRM 
combustion instability modeling, recasting the 
formulation into the traditional alpha notation. 
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In this formulation it is assumed that unsteady 
entropy is equal to zero. Both combustion stability 
models depend on an established quasi-steady and 
unsteady flow-field. These parameters are 
established using the COMSOL multiphysics 
finite element software. 

1.3 Unsteady Flow Calculations 

The above combustion stabilityframework do 
not attempt to model the entire time dependent 
flow-field of a SRM. This would require a costly 
and time consuming time accurate CFD model 
along with a yet to be defined acoustic admittance 
boundary conditions for the motor propellant. 
Traditionally the quasi-steady and unsteady flow- 
fields are solved separately, with the unsteady 
variables modeled via the homogenous wave 
equation. 
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This form of the wave equation assumes no 
mean flow contributions. These assumptions are 
valid for the majority of the rocket chamber, but 
break down within the converging section of the 
nozzle. The expulsion of unsteady energy through 
the nozzle of a rocket is identified as the 
predominate source of acoustic damping for most 
rocket systems. Recently, an approach to address 
nozzle damping with mean flow effects was 
implemented by French [7]. This new approach 
solves the acoustic velocity potential equation 
(AVPE) formulated by perturbing the Euler 
equations. 
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The AVPE includes effects of the mean flow 
variables on the unsteady acoustic mode shape 
and can be applied up to the sonic line. The 
present study utilizes both formulations as a way 
of verifying the COMSOL implementation 
against SSP and displaying the benefits of the 
more accurate AVPE. 

3. Use of COMSOL Multiphysics 

Implementing the above combustion 
instability models in COMSOL multiphysics 
requires the use of High Mach Number Flow 
(HMNF) CFD module, the Coefficient Form PDE 
module, and the Pressure Acoustics Module. First 
the HMNF module is used to model the SRM 
internal ballistics. The CFD results are post 


processed and the Mach equal to one surface is 
extracted. The Mach plane is used to generate a 
new geometry for the proceeding acoustic 
analyses, which do not model the nozzle exit 
cone. The results of the HMNF are used to define 
the mean flow terms in the AVPE. Results from 
the CFD and acoustic analyses are combined in 
post processing to calculate the stability alpha. 


3.1 High Mach Number Flow Analysis 


The present study employs COMSOL 
multiphysics finite element framework to model 
the steady flow-field parameters of a generic 
SRM using the Spalart-Allmaras turbulent flow 
equations within the HMNF module. A steady 
state analysis was performed with the 
initialization necessary for the Spalart-Allmara 
model. Equations 17-19 provide an overview of 
the conservation of mass, momentum and energy, 
respectively. 
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The HMNF module makes use of these fully 
compressible Navier- Stokes equation for an ideal 
gas, which are detailed in the COMSOL 
documentation. In order to model the injection of 
hot gas due to the burning propellant St. Robert’s 
law to model the solid propellant burn rate defined 


as, 

f = ap n (20) 

with f being the burn rate, p is the pressure at the 
wall, n is the burn rate exponent and a is an 
empirical factor. The normal velocity of the gas is 
then modeled as, 



with p p as the density of the propellant, and p g as 
the density of the combustion products. All other 
solid boundaries are modeled with the no- slip 
boundary condition, and the exit plane is modeled 
with the hybrid outflow condition. Figure 1 
displays the HMNF geometry with the flow 
boundary conditions identified. 
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Figure 1. SRM geometry with boundary conditions. 


3.2 Pressure Acoustics Analysis 

The Pressure Acoustics (PA) module is used 
to provide the unsteady flow parameters for 
combustion stability modeling. The PA module 
performs and Eigenvalue study on equation 14 
using the results from the HMNF analysis to 
supply the mean flow speed of sound and density. 
An acoustic hard wall boundary condition is used 
on all surfaces including the sonic (Mach = 1) 
line. Results from the Eigenvalue analysis and 
steady state analysis are combined in eq. (6) to 
calculate the motor’s growth alpha. The results 
are compared with the industry standard, SSP, as 
a verification of the COMSOL implementation. 
Figure 2 displays the geometry utilized in the 
unsteady Eigenvalue analyses. 

3.2 Coefficient Form Partial Differential 
Equation Analysis 

The Coefficient Form Partial Differential 
Equation module is used to determine the 
complex eigenvalues of the AVPE. Similarly to 
the PA analysis, the mean flow terms in eq. (16) 
are derived from the HMNF results. The 
mathematics model is truncated at the nozzle 
sonic line, where a zero flux boundary condition 
is self-satisfying. The remaining boundaries are 
modeled with a zero flux boundary condition, 
assuming zero acoustic absorption on all surfaces. 
The results of the steady-state CFD and AVPE 
analyses are used to calculate the linear acoustic 
growth rate as is defined by Flandro and Jacob [2, 

3]. 


4. Results 

The pertinent SRM ballistics parameters 
derived from the HMNF analysis are compared to 
the SPP results in table 2. The HMNF model 
shows a higher pressure at the head and aft end of 
the motor along with an elevated mass flow rate. 
These differences will have little effect on the 
stability analysis and comparison with SSP. It is 
expected that the differences arise due to how the 
initial guess used in the COMSOL analysis 
interacts with the inlet B.C. For instance, if the 
initial pressure was 8 14 psi for the entire chamber, 
then the burn rate, injection Mach number, and 
mass flow rate would be higher than if the initial 
guess had a linear pressure drop down the bore of 
the motor. But, the results indicate that, with some 
future work, the burning propellant of a SRM can 
be successfully modeled using the HMNF module 
in COMSOL. 

Table 2: Comparison of Ballistics Parameters 



Ph 

(psi) 

Pa 

(psi) 

771 

(Vo/s) 

Thrust 

(ib) 

HMNF 

834 

111 

98.84 

27130 

SPP 

818 

757 

95.15 

26690 

% diff 

1.95 

2.64 

3.88 

1.65 


A frequency comparison between the PA 
module and SSP is presented in table 3. Since the 
SSP analysis is strictly one dimensional, some 
differences are expected due to sudden area 
changes. Also, the geometry in the SSP model is 
slightly shorter than the COMSOL model due to 
its inability to model the sonic line. Figure 3 
displays the normalized acoustic mode shapes of 
the first six longitudinal modes derived by the PA 
analysis. 

Table 3: Comparison of Acoustic Frequencies 


Freq. 

(Hz) 

1L 

2L 

3L 

4L 

5L 

6L 

PA 

115 

231 

346 

462 

578 

695 

SSP 

116 

233 

350 

467 

584 

701 

% diff 

0.86 

0.86 

1.14 

1.07 

1.03 

0.86 



Figure 2. SRM geometry truncated at the sonic line. 
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Figure 3. Acoustic mode shapes derived from PA 
Eigenvalue analysis. 


6 displays the normalized imaginary part of the 
AVPE mode shapes. The imaginary portion of the 
mode shape represent a phase shift derived from 
mean flow effects. 


Alpha Comparison 
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Using the mode shapes from fig. 3 and the 
results from the HMNF analysis the right hand 
side of eq. (6) are evaluated and presented in fig. 
4. For each of the six longitudinal modes the 
growth rate is divided simplified into three terms; 
Pressure Coupling (PC), Flow Turning (FT), and 
Rotational Flow Correction (RFC). These terms 
represent the major contributors that can be 
compared in this analysis. Another major 
contributor that is absent is the Nozzle Damping 
(ND) term. The comparison would not be credible 
do to the differences between the two models. 
Figure 4 display a good comparison between the 
two approaches. The COMSOF results are 
consistently smaller than the SSP model, but do 
not diverge greatly. It should be noted that the FT 
alpha is actually a negative value, but is presented 
as a magnitude for ease of comparison. 

With a good comparison between two codes 
using the same stability model, focus is turned to 
a comparison of the homogeneous wave equation 
and AVPE. The inclusion of the mean flow terms 
on the right hand side of the wave equation act to 
shift the frequency and mode shape of the 
longitudinal acoustic modes. Figure 5 displays the 
normalized real component of the acoustic mode 
shapes for the first six longitudinal modes derived 
by the AVPE analysis. The results from the AVPE 
analysis have a real and imaginary component due 
to the acoustic/mean flow interactions. 

Comparing the mode shapes in figures 3 and 
5 little difference is observed in the first three 
acoustic modes, with noticeable variances in the 
latter three mode shapes. Focusing on the 
converging section of the nozzle on the right hand 
side of fig. 5 the effect of compressibility on the 
acoustic mode shapes is clearly observed. Figure 


Figure 4. Alpha comparison between SSP and 
COMSOL. 
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Figure 5. Real Acoustic mode shapes derived from 
AVPE Eigenvalue analysis. 
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Figure 6. Imaginary Acoustic mode shapes 
derived from AVPE Eigenvalue analysis. 


Along with a shift in the acoustic phasing, the 
inclusion of mean flow terms also shift the 


frequency of oscillation. A comparison of the PA 
frequencies and the AVPE frequencies is 
presented in table 4. Very little frequency shift is 
noticed in the low frequencies, with an increase in 
effect at higher frequencies. The subtle changes to 
the acoustic mode shape and frequency can have 
a large impact on the stability modeling. Acoustic 
results from PA and AVPE analyses are used to 
calculate PC and ND alpha terms using the energy 
based formulation. The bar graph in figure 7 
displays a comparison of the results. 


Table 4: Comparison of Acoustic Frequencies 


Freq. 

(Hz) 

1L 

2L 

3L 

4L 

5L 

6L 

PA 

115 

231 

346 

462 

578 

695 

AVPE 

115 

230 

345 

460 

576 

692 

% diff 

0.0 

0.43 

0.29 

0.43 

0.35 

0.43 


Alpha Comparison 
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■ PA-PC ■AVPE-PC ■ PA-ND ■ AVPE-ND 

Figure 7. Energy Corollary alpha comparison 
between PA and AVPE. 

5. Conclusions 

COMSOL multiphysics finite element 
software was used to model the combustion 
instability of a simplified SRM. The internal 
ballistics of the SRM were modeled with the 
HMNF CFD module and verified against an 
industry standard 1-D analysis tool. Two acoustic 
models were created using the Pressure Acoustics 
module and the Coefficient Form Partial 
Differential Equation module. The results from 
the PA module were used in conjunction with the 
CFD results to model the SRM stability using the 
Culick framework. These results were compared 
with the industry standard SSP for verification. 
Finally, the results from the AVPE were 
combined with the CFD results to model the SRM 
stability using the energy corollary framework. 
To demonstrate the advantage to using the AVPE, 


the stability alphas were calculated with the PA 
and math module solutions. 

This study has succeeded in establishing a 
proof of concept for using COMSOL 
multiphysics as a SRM combustion instability 
modeling tool. Future work will be focused on 
refining the inlet boundary condition modeled 
using the propellant bum rate. Also, further model 
will need to be developed to include vortical and 
entropy waves into the energy corollary. 
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